First-principles study of the phonon-limited mobility in n-type single-layer M0S2 
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In the present work we calculate the phonon-limited mobility in intrinsic n-type single-layer M0S2 
as a function of carrier density and temperature for T > 100 K. Using a first-principles approach 
for the calculation of the electron-phonon interaction, the deformation potentials and Frohlich in- 
teraction in the isolated M0S2 layer are determined. We find that the calculated room-temperature 



mobility of ~ 410 cm V 



is dominated by optical phonon scattering via deformation potential 



couplings and the Frohlich interaction with the deformation potentials to the intravalley homopo- 
lar and intervalley longitudinal optical phonons given by 4.1 x 10 s eV/cm and 2.6 x 10 s eV/cm, 
respectively. The mobility is weakly dependent on the carrier density and follows a /1 ~ T ' tem- 
perature dependence with 7 = 1.69 at room temperature. It is shown that a quenching of the 
characteristic homopolar mode which is likely to occur in top-gated samples, boosts the mobility 
with ~ 70 cm V -1 s" 1 and can be observed as a decrease in the exponent to 7 = 1.52. Our 
findings indicate that the intrinsic phonon-limited mobility is approached in samples where a high-K 
dielectric that effectively screens charge impurities is used as gate oxide. 
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PACS numbers: 81.05.Hd, 72.10.-d, 72.20.-i, 72.80.Jc 



I. INTRODUCTION 



Alongside with the rise of graphene and the explo- 
ration of its unique electronic properties^— the search 
for other two-dimensional (2D) materials with promis- 
ing electronic properties has gained increased interest. 4 
Metal dichalcogenides which are layered materials similar 
to graphite, provide interesting candidates. Their layered 
structure with weak inter/layer van der Waals bonds al- 
lows for fabrication of single- to few-layer samples using 
mechanical peeling/cleavage or chemical exfoliation tech- 
niques similar to the fabrication of graphene^ 8 How- 
ever, in contrast to graphene they are semiconductors 
and hence come with a naturally occurring band gap — a 
property essential for electronic applications. 

The vibrational and optical properties of single- to few- 
layer samples of M0S2 have recently been studied ex- 
tensively using various experimental techniques^— In 
contrast to the bulk material which has an indirect-gap, 
it has been demonstrated that single-layer M0S2 is a 
direct-gap semiconductor with a gap of 1.8 eV. 10 ' 12 To- 
gether with the excellent electrostatic control inherent 
of two-dimensional materials, the large band gap makes 
it well-suited for low power applications.— So far, elec- 
trical characterizations of single-layer M0S2 have shown 
n-type conductivity with room-temperature mobilities in 
the range 0.5 — 3 cm 2 V -1 s -1 ^^ Compared to early 
studies of the intralayer mobility of bulk M0S2 where 
mobilities in the range 100 — 260 cm 2 V -1 s _1 were 
reported^ this is rather low. In a recent experiment 
the use of a high-K gate dielectric in a top-gated device 
was shown to boost the carrier mobility to a value of 
200 cm 2 V -1 s _1 4 The observed increase in the mobility 
was attributed to screening of impurities by the high-K 
dielectric and/or modifications of M0S2 phonons in the 
top-gated sample. 



In order to shed further light on the measured mobili- 
ties and the possible role of phonon damping due to a top- 
gate dielectric^ we have carried out a detailed study of 
the phonon-limited mobility in single-layer MoS2- Since 
phonon scattering is an intrinsic scattering mechanism 
that often dominates at room temperature, the theo- 
retically predicted phonon-limited mobility sets an up- 
per limit for the experimentally achievable mobilities in 
single-layer M0S2. 15 In experimental situations, the tem- 
perature dependence of the mobility can be useful for the 
identification of dominating scattering mechanisms and 
also help to establish the value of the intrinsic electron- 
phonon couplings. 16 However, the existence of additional 
extrinsic scattering mechanisms often complicates the in- 
terpretation. For example, in graphene samples impu- 
rity scattering and scattering on surface polar optical 
phonons of the gate oxide are important mobility-limiting 
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FIG. 1: (Color online) Atomic structure and conduction band 
of single-layer M0S2. Left: Primitive unit cell and structure 
of an M0S2 layer with the molybdenum and sulfur atoms po- 
sitioned in displaced hexagonal layers. Right: Contour plot 
showing the lowest lying conduction band valleys as obtained 
with DFT-LDA in the hexagonal Brillouin zone of single-layer 
M0S2. For n-type M0S2, the low-field mobility is determined 
by the carrier properties in the K, Tf'-valleys which are sepa- 
rated in energy from the satellite valleys inside the Brillouin 
zone. 
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factorsji^— It can therefore be difficult to establish the 
value and nature of the intrinsic electron-phonon interac- 
tion from experiments alone. Theoretical studies there- 
fore provide useful information for the interpretation of 
experimentally measured mobilities. 

In the present work, the electron-phonon interaction in 
single-layer M0S2 is calculated from first-principles using 
a density-functional based approach. From the calcu- 
lated electron-phonon couplings, the acoustic (ADP) and 
optical deformation potentials (ODP) and the Frohlich 
interaction in single-layer M0S2 are inferred. Using 
these as inputs in the Boltzmann equation, the phonon- 
limited mobility is calculated in the high-temperature 
regime (T > 100 K) where phonon scattering often dom- 
inate the mobility. Our findings demonstrate that the 
calculated phonon-limited room-temperature mobility of 
~ 410 cm 2 V -1 s _1 is dominated by optical deformation 
potential and polar optical scattering via the Frohlich 
interaction. The dominating deformation potentials of 
L>°p = 4.1 x 10 8 eV/cm and Z)£ = 2.6 x 10 8 eV/cm 
originate from the couplings to the homopolar and the 
intervalley polar LO phonon. This is in contrast to pre- 
vious studies for bulk M0S2, where the mobility has been 
assumed to be dominated entirely by scattering on the 
homopolar mode . 14 ' 23 ! 24 Furthermore, we show that a 
quenching of the characteristic homopolar mode which is 
polarized in the direction normal to the M0S2 layer can 
be observed as a change in the exponent 7 of the generic 
temperature dependence fi ~ 7 1-7 of the mobility. Such a 
quenching can be expected to occur in top-gated samples 
where the M0S2 layer is sandwiched between a substrate 
and a gate oxide. 

The paper is organized as follows. In Section Ql] the 
band structure and phonon dispersion of single-layer 
M0S2 as obtained with DFT-LDA are presented. The lat- 
ter is required for the calculation of the electron-phonon 
coupling presented in Section IIIII From the calculated 
electron-phonon couplings, we extract zero- and first- 
order deformation potentials for both intra and inter- 
valley phonons. In addition, the coupling constant for 
the Frohlich interaction with the polar LO phonon is de- 
termined. In Section IIV1 the calculation of the phonon- 
limited mobility within the Boltzmann equation is out- 
lined. This includes a detailed treatment of the phonon 
collision integral from which the scattering rates for the 
different coupling mechanisms can be extracted. In Sec- 
tion [V] we present the calculated mobilities as a function 
of carrier density and temperature. Finally, in Section lVll 
we summarize and discuss our findings. 



II. ELECTRONIC STRUCTURE AND PHONON 
DISPERSION 

The electronic structure and phonon dispersion of 
single-layer M0S2 have been calculated within DFT in 
the LDA approximation using a real-space projector- 
augmented wave (PAW) method^ - — The equilibrium 
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FIG. 2: (Color online) Band structure of single-layer M0S2. 
The inset shows a zoom of the conduction band along the path 
T-K-M. The parabolic band (red dashed line) with effective 
mass m* = 0.48 m e is seen to give a perfect description of 
the conduction band valley in the if-point for the range of 
energies relevant for the low-field mobility. 

lattice constant of the hexagonal unit cell shown in Fig. Q] 
was found to be a = 3.14 A using a 11 x 11 k-point 
sampling of the Brillouin zone. In order to eliminate in- 
teractions between the layers due to periodic boundary 
conditions, a large interlayer distance of 10 A has been 
used in the calculations. 



A. Band structure 

In the present study, the band structure of single-layer 
MoS 2 is calculated within DFT-LDA, 2 ^ While DFT- 
LDA in general underestimates band gaps, the resulting 
dispersion of individual bands, i.e. effective masses and 
energy differences between valleys, is less problematic. 
We note, however, that recent GW quasi-particle cal- 
culations which in general give a better description of 
the band structure in standard semiconductors^ sug- 
gest that the distance between the K, if '-valleys and the 
satellite valleys located on the T-K path inside the Bril- 
louin zone not as large as predicted by DFT (see below) 
and that the valley ordering may even be invertedi 53 i 54 
This, however, seems to contradict with the experimen- 
tal consensus that single-layer M0S2 is a direct-gap semi- 
conducto r 10 ' 12 and further studies are needed to clarify 
this issue. In the following we therefore use the DFT- 
LDA band structure and comment on the possible con- 
sequences of a valley inversion in Section IVI1 

The calculated DFT-LDA band structure is shown in 
Fig. [2j It predicts a direct band gap of 1.8 eV at the 
if-poing in the Brillouin zone. The inset shows a zoom 
of the bottom of the conduction band along the T-K- 
M path in the Brillouin zone. As demonstrated by the 
dashed line, the conduction band is perfectly parabolic in 
the if- valley. Furthermore, the satellite valley positioned 
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on the path between the T-point and the if -point lies on 
the order of ~ 200 meV above the conduction band edge 
and is therefore not relevant for the low-field transport. 

The right plot in Fig.[T]shows a contour plot of the low- 
est lying conduction band valleys in the two-dimensional 
Brillouin zone. Here, the K, if' -valleys that are popu- 
lated in n-type M0S2, are positioned at the corners of 
the hexagonal Brillouin zone. Due to their isotropic and 
parabolic nature, the part of the conduction band rele- 
vant for the low-field mobility can to a good approxima- 
tion be described by simple parabolic bands 



with an effective electron mass of m* — 0.48 m e and 
where k is measured with respect to the K, if' -points 
in the Brillouin zone. The two-dimensional nature of the 
carriers is reflected in the constant density of states given 
by Po = g s g v m* where g s = 2 and g v — 2 are the 

spin and K, if '-valley degeneracy, respectively. The large 
density of states which follows from the high effective 
mass of the conduction band in the K, if '-valleys, results 
in non-degenerate carrier distributions except for very 
high carrier densities. 

B. Phonon dispersion 

The phonon dispersion has been obtained with the 
supercell-based small-displacement method 3 ^ using a 9 x 
9 supercell. The resulting phonon dispersion shown in 
Fig. [3] is in excellent agreement with recent calculations 
of the lattice dynamics in two-dimensional MoS2^ 

With three atoms in the unit cell, single-layer M0S2 
has nine phonon branches — three acoustic and six opti- 
cal branches. Of the three acoustic branches, the fre- 
quency of the out-of-plane flexural mode is quadratic in 
q for q — > 0. In the long-wavelength limit, the frequency 
of the remaining transverse acoustic (TA) and longitudi- 
nal acoustic (LA) modes are given by the in-plane sound 
velocity c\ as 

WqA = c\q. (2) 

Here, the sound velocity is found to be 4.4 x 10 3 and 
6.5 x 10 3 m/s for the TA and LA mode, respectively. 

The gap in the phonon dispersion completely sepa- 
rates the acoustic and optical branches even at the high- 
symmetry points at the zone-boundary where the acous- 
tic and optical modes become similar. The two lowest 
optical branches belong to the non-polar optical modes. 
Due to an insignificant coupling to the charge carriers, 
they are not relevant for the present study. 

The next two branches with a phonon energy of ~ 
48 meV at the T-point are the transverse (TO) and lon- 
gitudinal (LO) polar optical modes where the Mo and S 
atoms vibrate in counterphase. In bulk polar materials, 
the coupling of the lattice to the macroscopic polariza- 
tion setup by the lattice vibration of the polar LO mode 
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FIG. 3: Phonon dispersion of single-layer M0S2 calculated 
with the small displacement method (see e.g. Ref.30) using 
a 9 x 9 supercell. The frequencies of the two optical Raman 
active E2 S and Ai 9 modes at 48 meV and 50 meV, respec- 
tively, are in excellent agreement with recent experimental 
measurements.— 



results in the so-called LO-TO splitting between the two 
modes in the long- wavelength limit. The inclusion of 
this effect from first-principles requires knowledge of the 
Born effective charges.— ~— In two-dimensional materials, 
however, the lack of periodicity in the direction perpen- 
dicular to the layer removes the LO-TO splitting. 35 The 
coupling to the macroscopic polarization of the LO mode 
will therefore be neglected here. 

The almost dispersionless phonon at ~ 50 meV is the 
so-called homopolar mode which is characteristic for lay- 
ered structures. The lattice vibration of this mode cor- 
responds to a change in the layer thickness and has the 
sulfur layers vibrating in counterphase in the direction 
normal to the layer plane while the Mo layer remains 
stationary. The change in the potential associated with 
this lattice vibration has previously been demonstrated 
to result in a large deformation potential in bulk MoS2- 14 



III. ELECTRON-PHONON COUPLING 

In this section, we present the electron-phonon cou- 
pling in single-layer M0S2 obtained with a first-principles 
based DFT approach. The method is based on a su- 
percell approach analogous to that used for the calcula- 
tion of the phonon dispersion and is outlined in App. [X] 
The calculated electron-phonon couplings are discussed 
in the context of deformation potential couplings and the 
Frohlich interaction well-known from the semiconductor 
literature. 

Within the adiabatic approximation for the electron- 
phonon interaction, the coupling strength for the phonon 
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mode with wave vector q and branch index A is given by 



9 ^ V 2MN^ M ^ (3) 

where cj q A is the phonon frequency, M is an appropriately 
defined effective mass, N is the number of unit cells in 
the crystal, and 

M^ q =(k + q|(5V qA (r)|k) (4) 

is the coupling matrix element where k is the wave vector 
of the carrier being scattered and 5Vq\ is the change in 
the effective potential per unit displacement along the 
vibrational normal mode. 

Due to the valley degeneracy in the conduction band, 
both intra and intervalley phonon scattering of the carri- 
ers in the K, if'-valleys need to be considered. Here, the 
coupling constants for these scattering processes are ap- 
proximated by the electron-phonon coupling at the bot- 
tom of the valleys, i.e. with k = K,K'. With this ap- 
proach, the intra and intervalley scattering for the K, K'- 
valleys are thus assumed independent on the wave vector 
of the carriers. 

In the following two sections, the calculated electron- 
phonon couplings are presented^ The different deforma- 
tion potential couplings are discussed and the functional 
form of the Frohlich interaction in 2D materials is es- 
tablished. Piezoelectric coupling to the acoustic phonons 
which occurs in materials without inversion symmetry, 
is most important at low temperatures^ and will not 
be considered here. As deformation potentials are of- 
ten extracted as empirical parameters from experimental 
mobilities, it can be difficult to disentangle contributions 
from different phonons. Here, the first-principles calcula- 
tion of the electron-phonon coupling allows for a detailed 
analysis of the couplings and to assign deformation po- 
tentials to the individual intra and intervalley phonons. 

A. Deformation potentials 

The deformation potential interaction describes how 
carriers interact with the local changes in the crystal po- 
tential associated with a lattice vibration. Within the de- 
formation potential approximation, the electron-phonon 
coupling is expressed as2& 

where A is the area of the sample, p is the atomic mass 
density per area and M q A is the coupling matrix element 
for a given valley which is assumed independent on the 
k- vector of the carriers. This expression follows from 
the general definition of the electron-phonon coupling in 
Eq. [3] by setting the effective mass M equal to the sum of 
the atomic masses in the unit cell. With this convention 



for the effective mass, MN — Ap, and the expression in 
Eq. ([5]) is obtained. 

For scattering on acoustic phonons, the coupling ma- 
trix element is linear in q in the long- wavelength limit, 

M q A = S A <7, (6) 

where Sa is the acoustic deformation potential. 

In the case of optical phonon scattering, both coupling 
to zero- and first-order in q must be considered. The in- 
teraction via the constant zero-order optical deformation 
potential D x is given by 

MqA = D° x . (7) 

The coupling via the zero-order deformation potential is 
dictated by selection rules for the coupling matrix ele- 
ments. Therefore, only symmetry-allowed phonons can 
couple to the carriers via the zero-order interaction. The 
coupling via the first-order interaction is given by Eq. ([5]) 
with the acoustic deformation potential replaced by the 
first-order optical deformation potential D\. Both the 
zero- and first-order deformation potential coupling can 
give rise to intra and intervalley scattering. 

The absolute value of the calculated coupling matrix el- 
ements |M q A| are shown in Fig. U for phonons which cou- 
ple to the carriers via deformation potential interactions. 
They are plotted in the range of phonon wave vectors rel- 
evant for intra and intervalley scattering processes^ and 
only phonon modes with significant coupling strengths 
have been included. Although the coupling matrix ele- 
ments are shown only for k = K, the matrix elements 
for k = K' are related through time-reversal symmetry 
as M^ x = M_„ x . The three-fold rotational symmetry 
of the coupling matrix elements in q-space stems from 
the symmetry of the conduction band in the vicinity of 
the K, if '-valleys (see Fig.[TJ|. Since the symmetry of the 
matrix elements has not been imposed by hand, slight 
deviations from the three-fold rotational symmetry can 
be observed. 

The acoustic deformation potential couplings for the 
TA and LA modes are shown in the two top plots of 
Fig. |4j Due to the inclusion of Umklapp processes here, 
the coupling to the TA mode does not vanish as is most 
often assumed4£ Only along high-symmetry directions in 
the Brillouin zone is this the case. This results in a highly 
anisotropic coupling to the TA mode. On the other hand, 
the deformation potential coupling for the LA mode is 
perfectly isotropic in the long-wavelength limit but also 
becomes anisotropic at shorter wavelengths. In agree- 
ment with Eq. (j6)), both the TA and LA coupling matrix 
elements are linear in q in the long- wavelength limit. 

The next two plots show the couplings to the intraval- 
ley polar TO and homopolar optical modes. While the 
interaction with the polar TO phonon corresponds to 
a first-order optical deformation potential, the interac- 
tion with the homopolar mode acquires a finite value of 
~ 4.8 eV/A in the F-point and corresponds to a strong 
zero-order deformation potential coupling. The large 
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FIG. 4: (Color online) Deformation potential couplings in 
single-layer M0S2. Only phonon modes with significant cou- 
pling strength are shown. The contour plots show the absolute 
value of the coupling matrix elements |M q x| in Eq. Q with 
k = K and k' = K + q as a function of the two-dimensional 
phonon wave vector q (note the different color scales in the 
plots). The upper four plots correspond to intravalley scat- 
tering and the two bottom plots to intervalley scattering on 
phonons with wave vector q = 2K. The two zero-order defor- 
mation potential couplings to the intravalley homopolar and 
intervalley LO phonon play an important role for the room- 
temperature mobility. 



have significant weight here, the homopolar lattice vi- 
bration gives rise to a significant shift of the K, if' -valley 
states. The large deformation potential associated with 
the coupling to the homopolar mode is also present in 
bulk MoS 2 ^ 

The two bottom plots in Fig. @] show the couplings 
for the intervalley TA and polar LO phonons. Due to 
the nearly constant frequency of the intervalley acoustic 
phonons, their couplings are classified as optical defor- 
mation potentials in the following. Both the shown cou- 
pling to the intervalley TA phonon and the coupling to 
the intervalley LA phonon result in first-order deforma- 
tion potentials. The intervalley coupling for the polar 
LO phonon gives rise to a zero-order deformation poten- 
tial which results from an identical in-plane motion of the 
two sulfur layers in the lattice vibration of the intervalley 
phonon. 

In general, the deformation potential approximation, 
i.e. the assumption of isotropic and constant/linear cou- 
pling matrix elements, in Eqs. © and (JTJ) is seen to 
hold in the long-wavelength limit only. At shorter wave- 
lengths, the first-principles couplings become anisotropic 
and have a more complicated q-dependence. When de- 
termined experimentally from e.g. the temperature de- 
pendence of the mobility^ the deformation potentials 
in Eqs. ([5]) and J7]) implicitly account for the more com- 
plex q-dependence of the true coupling matrix element. 
In Section |Vj the theoretical deformation potentials for 
single-layer M0S2 are determined from the first-principles 
electron-phonon couplings. This is done by fitting their 
associated scattering rates to the scattering rates ob- 
tained with the first-principles coupling matrix elements. 
The resulting the deformation potentials, can be used in 
practical transport calculations based on e.g. the Boltz- 
mann equation or Monte Carlo simulations. It should 
be noted that similar routes for first-principles calcula- 
tions of deformation potentials have been given in the 
literature4ir— 



B. Frohlich interaction 

The lattice vibration of the polar LO phonon gives 
rise to a macroscopic electric field that couples to the 
charge carriers. For bulk three-dimensional systems, the 
coupling to the field is given by the Frohlich interaction 
which diverges as 1/q in the long- wavelength limit^ 



deformation potential coupling to the homopolar mode 
stems from its characteristic lattice vibration polarized 
in the direction perpendicular to the layer corresponding 
to a change in the layer thickness. The associated change 
in the effective potential from the counterphase oscilla- 
tion of the two negatively charged sulfur layers results 
in a significant change of the potential towards the cen- 
ter of the M0S2 layer. As the electronic Bloch functions 



1 e^hoj ( 1 1 \ 

9w{q) = ~ q {-^ {-^ - V^) ' (8) 

where e is the vacuum permittivity, V is the volume of 
the sample , and e°° and e° are the high-frequency optical 
and static dielectric constant, respectively. 

In atomically thin materials, the two-dimensional na- 
ture of both the LO phonon and the charge carriers leads 
to a qualitatively different g-dependence of the Frohlich 
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FIG. 5: (Color online) Frohlich interaction for the polar opti- 
cal LO mode. The dashed line shows the analytical coupling 
in Eq. ((9| with the coupling constant gFi = 98 meV and 
the effective layer thickness a = 4.41 A fitted to the long- 
wavelength limit of the calculated coupling (dots). 



interaction. The situation is similar to 2D semiconduc- 
tor heterostructures where the Frohlich interaction has 
been studied using dielectric continuum models and mi- 
croscopic based approaches From the microscopic 
considerations presented in App. [B] we derive the fol- 
lowing functional form for the Frohlich interaction to the 
polar LO phonon in 2D materials, 



3lo(<?) = ffFr x erfc(gcr/2). 



(9) 



Here, the coupling constant is the equivalent of the 
square root factors in Eq. ©, a is the effective width of 
the electronic Bloch states and erfc is the complementary 
error function. 

Figure [5] shows the first-principles coupling <7lo to the 
polar LO phonon in single-layer M0S2 (dots) as a func- 
tion of the phonon wave vector along the T-K path. The 
dashed line shows a fit of the coupling in Eq. to the 
long-wavelength limit of the calculated coupling. With a 
coupling constant of <?Fr = 98 meV and an effective width 
of a = 4.41 A, the analytic form gives a perfect descrip- 
tion of the calculated coupling in the long-wavelength 
limit. For shorter wavelengths, the zero-order deforma- 
tion potential coupling to the intervalley LO phonon from 
Fig. [4] becomes dominant, and a deviation away from the 
behavior in Eq. ^ is observed. 



IV. BOLTZMANN EQUATION 

In the regime of diffuse transport dominated by phonon 
scattering, the mobility can be obtained with semiclas- 
sical Boltzmann transport theory. In the absence of 
spatial gradients the Boltzmann equation for the out-of- 
equilibrium distribution function / k of the charge carriers 
reads^ 



dfu ■ <9/ k 
^ + k.V k / k = — 



(10) 



where the time evolution of the crystal momentum is gov- 
erned by Tik = q~E and q is the charge of the carriers. In 
the case of a time-independent uniform electric field, the 
steady-state version of the linearized Boltzmann equation 
takes the form 



<l 



E 



V k / k ° 



k B T~ 
dt 



E-v k /°(l-/°) 



(11) 



where = / (e k ) is the equilibrium Fermi-Dirac dis- 
tribution function, df£/de k = -/^(l - f£)/k B T, and 
v k = V k e k /fi. is the band velocity. In linear response, 
the out-of-equilibrium distribution function can be writ- 
ten as the equilibrium function plus a small deviation Sf 
away from its equilibrium valued i.e. 



5h = /k - /£ = / k (i - £)* 



(12) 



Here, the last equality has denned the deviation func- 
tion V'k- Furthermore, following the spirit of the itera- 
tive method of Rode^ the angular dependence of the 
deviation function is separated out as 



ipu = <\>k cos ( 



(13) 



where 6> k is the angle between the k-vector and the force 
exerted on the carriers by the applied field E. In gen- 
eral, (f>k is a function of the k-vector and not only its 
magnitude k. However, for isotropic bands as considered 
here, the angular dependence of the deviation function is 
entirely accounted for by the cosine factor in Eq. (fT5|) . 38 
Considering the phonon scattering, the collision inte- 
gral describes how accelerated carriers are driven back 
towards their equilibrium distribution by scattering on 
acoustic and optical phonons. In the high-temperature 
regime of interest here, the collision integral can be split 
up in two contributions. The first, accounting for the 
quasielastic scattering on acoustic phonons, can be ex- 
pressed in the form of a relaxation time r e i. The second, 
describing the inelastic scattering on optical phonons, 
will in general be an integral operator 2i nc i. The colli- 
sion integral can thus be written 



dt 



coll 



^ k , T 111 



(14) 



where the explicit forms of the relaxation time and the 
integral operator will be considered in detail below. 

With the above form of the collision integral, the Boltz- 
mann equation can be solved by iterating the equation 



eEvk 
k B T 



1 



(15) 



for the deviation function Here, the angular depen- 
dence cos0 k and a /^(l — /£) factor have been divided 
out and the new integral operator I' incl is defined by 



-^incl 



/k°(l 



/k)cOS0k4el- 



(16) 
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From the solution of the Boltzmann equation the cur- 
rent and drift mobility of the carriers can be obtained. 
Taking the electric field to be oriented along the re- 
direction, the current density is given by 



v xx E x = ~~^E v xSfk, 

k 



(17) 



where a xx is the conductivity, A is the area of the sample, 
and Vi = hki /m* is the band velocity for parabolic bands 
with i = x,y, z. From the definition [i xx = <j xx /ne of the 
drift mobility, it follows that 



/"a 



where the modified deviation function 
time is defined by 



(18) 

t>k having units of 



- _ m*k B T 
0fe ~~ ~eEJik 

and the energy- weighted average (•) is defined by 



(A) 



1 



sfck p{eu.)e k A k - 



df 
<9e k 



(19) 



(20) 



Here, n is the two-dimensional carrier density. In the re- 
laxation time approximation 4>k — ifc and the well-known 
Drude expression for the mobility \i = e{Tk)/m* is recov- 
ered. 



A. Phonon collision integral 



The phonon collision integral has been considered in 
great detail in the literature for two-dimensional elec- 
tron gases in semiconductor heterostructure (see e.g. 
Ref|37J). In general, these treatments consider scattering 
on three-dimensional or quasi two-dimensional phonons. 
In atomically thin materials the phonons are strictly two- 
dimensional which results in a slightly different treat- 
ment. In the following, a full account of the two- 
dimensional phonon collision integral is given for scat- 
tering on acoustic and optical phonons. 

With the distribution function written on the form 
in Eq. (fl2|. the linearized collision integral for electron- 
phonon scattering takes the forrrt^ 



dt 



2tt 



IffqAl 2 



E 



/k (1 - /k+q) ^qA (tfk " ^k+q) <S(e k+q - 6k - 

f /k (1 - /k-J (1 + ^a) (</>k " </>k-q) ^k-q - £k + ^qA) 



(21) 



where N® x = N (hujq\) is the equilibrium distribution 
of the phonons given by the Bose-Einstein distribution 
function and /k± q = /°(£k ± fejq) is understood. The 
different terms inside the square brackets account for 
scattering out of (V'k) and into (V'kiq) the state with 
wave vector k via absorption and emission of phonons. 
Screening of the electron-phonon interaction by the car- 
riers themselves is accounted for by the static dielectric 
function e, which is both wave vector q and tempera- 
ture T dependent. Due to the large effective mass of 
the conduction band, the charge carriers in single-layer 
M0S2 are non-degenerate except at very high carrier con- 
centrations. Semiclassical screening where the screening 
length is given by the inverse of the Debye-Hiickel wave 
vector fcn = e 2 n/2k-QT therefore applies4£ For the con- 
sidered values of carrier densities and temperatures, this 
corresponds to a small fraction of the size of the Bril- 
louin zone, i.e. fcn -C 2n/a. As scattering on phonons 
in general involves larger wave vectors, screening by the 
carriers can to a good approximation be neglected here. 
The calculated mobilities therefore provide a lower limit 



for the proper screened mobility. 



1. Quasielastic scattering on acoustic phonons 

For quasielastic scattering on acoustic phonons at high 
temperatures, the energy of the acoustic phonon can be 
neglected in the collision integral in Eq. (|2"Tj) . As a con- 
sequence, the collision integral can be recast in the form 
of the following relaxation time^i 



E 

k' 



1 - cos6> kk /)P kk / 



(22) 



where the summation variable has been changed to k' = 
k ± q and the transition matrix element for quasielastic 
scattering on acoustic phonons is given by 

W = jEW 2 \N° x 5(e k ,-e k ) 



Pi 



;i + iv q o A ) ( 5( ek ,- £k ) 



(23) 
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For isotropic scattering as assumed here in Eq. ([6]), the 
square of the electron-phonon coupling can be expressed 



2Apc x 



(24) 



where the acoustic phonon frequency has been expressed 
in terms of the sound velocity c\ . Except at very low tem- 
peratures Sw q fceT implying that the equipartition ap- 
proximation N® ~ kBT/hujq 1 for the Bose- Einstein 
distribution applies. With the resulting q-factors in 
the transition matrix element in Eq. (|23j) canceling, the 
cos $kk' in the k'-sum in Eq. (|2"2"|) vanishes and the first 
term yields a factor density of states divided by the spin 
and valley degeneracies po / g s 9v ■ The relaxation time for 
acoustic phonon scattering becomes 



1 _ m*E 2 k B T 



W>pc\ 



(25) 



The independence on the carrier energy and the r^ 1 ~ T 
temperature dependence of the acoustic scattering rate 
are characteristic for charge carriers in two-dimensional 
heterostructures and layered materials ! 14 ' 37 



2. Inelastic scattering on dispersionless optical phonons 



scattering is given by 

Pkk ' = tE'^i 2 x 



N%6(ev - E k - hoj x ) 

+ (1 + AT°) 5(e k , - e k + Hlux) 

(28) 



For the in-scattering part, the desired cos 0k factor in 
Eq. (fTo]) can be extracted from -0k' using the rela- 
tion cos 0k' = cos 0k cos 0kk' — sin 0k sin 0kk' ■ Since the 
sine term vanishes from symmetry consideration, the in- 
scattering part of the inelastic collision integral reduces 
to 



k' 



COS k k' J Pkk' 



l-/k° 



(29) 



The evaluation of the k'-sum in Eqs. (|27p and (|2T)1) is 
outlined in App. [Cj Here, the assumption of dispersion- 
less optical phonons allows for an semianalytical treat- 
ment. For zero- and first-order coupling within the defor- 
mation potential approximation, the sum can be carried 
out analytically. The resulting expressions for the colli- 
sion i ntegr al are given in Eqs. (IC14I) . (|C15|) . (|C16|) . (|C17|) 
and (IC18[) . In the case of the Frohlich interaction, the 
angular part of the q-integral must be done numerically. 



For inelastic scattering on optical phonons the phonon 
energy can no longer be neglected and the collision in- 
tegral in Eq. ([21]) must be considered in full detail. Un- 
der the reasonable assumption of dispersionless optical 
phonons lu^x = uj\, the inelastic collision integral can be 
treated semianalytically. The overall procedure for the 
evaluation of the collision integral is given below, while 
the calculational details have been collected in App. [C] 
The resulting expressions for the collision integral apply 
to both intravalley and intervalley optical and intervalley 
acoustic phonons. 

In the following, the integral operator for inelastic scat- 
tering in Eq. (1161) is split up in separate out- and in- 
scattering contributions, 

i\nM = iZoAM + iZJM, (26) 

which include the terms in Eq. (|21[) involving ip^ and 
-0k±q, respectively. With the contributions from differ- 
ent phonon branches A adding up, scattering on a single 
phonon with branch index A is considered in the follow- 
ing. 

From Eq. ([2"Tj) . the out-scattering part of the collision 
integral follows directly as 

k' 1 J k 

where the transition matrix element for optical phonon 



3. Optical deformation potential scattering rates 

In spite of the fact that the collision integral for in- 
elastic scattering on optical phonons cannot be recast in 
the form of a (momentum) relaxation time^ a scatter- 
ing rate related to the inverse carrier lifetime can still 
be defined from the out-scattering part of the inelastic 
collision integral alone. The scattering rate so defined is 
given by 

i = E p ^'T3# ( 3 °) 

T k k ' 1 J k 

and corresponds to the imaginary part of the elec- 
tronic self-energy in the Born approximation. 46 Below, 
the resulting scattering rates for zero-order and first- 
order deformation potential scattering are given for non- 
degenerate carriers, i.e. with the Fermi factors in Eq. (j3TJ)) 
neglected. They follow straight-forwardly from the ex- 
pressions for the out-scattering part of the collision in- 
tegral derived in App. [C] It should be noted that the 
scattering rate for the zero-order deformation poten- 
tial interaction given below in Eq. (|3f j) . in fact defines 
a proper momentum relaxation time because the in- 
scattering part of the collision integral vanishes in this 
case (see App. ICj). 

For zero-order deformation potential scattering, the 
scattering rate is independent of the carrier energy and 
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given by 



1 



m*(Dlf 



^incl 



1 + e ritJA / feBT e(e k - huj x ) 



N° x . (31) 



Here, Q(x) is the Heavyside step function which assures 
that only electrons with sufficient energy can emit a 
phonon. 

The scattering rate for coupling via the first-order de- 
formation potential is found to be 



1 



,-inel 



m* 2 {D\f 



(2e k + Tvljx) + e h ^ k * T 
x 6(e k - hhj x ) (2e k - ftw x ) 



(32) 



Due to the linear dependence on the carrier energy, zero- 
order scattering processes dominate first-order processes 
at low energies. Only under high-field conditions where 
the carriers are accelerated to high velocities will first- 
order scattering become significant. 

The expressions for the scattering rates in Eqs. (|3Tj) 
and (|32[) above apply to scattering on dispersionless in- 
tervalley acoustic phonons and intra/intervalley optical 
phonons. Except for a factor of \/ e-^ ± hu\ originat- 
ing from the density of states, the energy dependence 
of the scattering rates is identical to that of their three- 
dimensional analogs^ 



V. RESULTS 

In the following, the scattering rate and pho non- 
limited mobility in single-layer M0S2 are studied as a 
function of carrier energy, temperature T and carrier den- 
sity n using the material parameters collected in Tab. [I] 
Here, the reported deformation potentials represent ef- 
fective coupling parameters for the deformation potential 
approximation in Eqs. (0 and (see below). 



A. Scattering rates 

With access to the first-principles electron-phonon 
couplings, the scattering rate taking into account the 
anisotropy and more complex g-dependence of the first- 
principles coupling matrix elements can be evaluated. 
They are obtained using the expression for the scattering 
rate in Eq. (|30[) with the respective transition matrix el- 
ements for acoustic and optical phonon scattering given 
in Eqs. (|2"3")l and (|2"51). 50 In order to account for the dif- 
ferent coupling matrix elements in the K, K'- valleys, the 
scattering rates have been averaged over different high- 
symmetry directions of the carrier wave vector k. The 
resulting scattering rates are shown in Fig. [5] as a func- 
tion of the carrier energy for non-degenerate carriers at 
T = 300 K. The different lines show the contributions to 



the total scattering rate from the various electron-phonon 
couplings to the intra and intervalley phonons which have 
been grouped according to their coupling type, i.e. acous- 
tic deformation potentials (ADPs), zero/first-order opti- 
cal deformation potentials (ODPs) and the Frohlich in- 
teraction. The acoustic deformation potential scatter- 
ing includes the quasielastic intravalley scattering on the 
TA and LA phonons with linear dispersions. Scatter- 
ing on intervalley acoustic phonons is considered as op- 
tical deformation potential scattering. Both the total 
scattering rate and the contributions from the different 
coupling types have been obtained using Matthiessen's 
rule by summing the scattering rates from the individual 
phonons, i.e. 



St 1 - 



(33) 



For carrier energies below the optical phonon frequen- 
cies, the total scattering rate is dominated by acoustic de- 
formation potential scattering. At higher energies zero- 
order deformation potential scattering and polar optical 
scattering via the Frohlich interaction become dominant. 
Due to the linear dependence on the carrier energy, the 
first-order deformation potential scattering on the inter- 
valley acoustic phonons and the optical phonons is in 



Parameter 



Symbol 



Value 



Lattice constant 
Ion mass density 
Effective electron mass 
Transverse sound velocity 
Longitudinal sound velocity 
Acoustic deformation potentials 
TA 
LA 

Optical deformation potentials 

TA 

LA 

TO 

TO 

LO 

Homopolar 

Frohlich interaction (LO) 
Effective layer thickness 
Coupling constant 
Optical phonon energies 
Polar LO 
Homopolar 



a 
P 

m* 

CTA 

cla 

Hta 
Hla 



^K,TA 
£*K,LA 
^r,TO 
fk,TO 
^K.LO 



a 
/kj H p 



3.14 A 
3.1 x 1(T 7 g/cm 2 
0.48 m e 
4.2 x 10 3 m/s 
6.7 x 10 3 m/s 

1.6 eV 

2.8 eV 

5.9 eV 
3.9 eV 
4.0 eV 
1.9 eV 

2.6 x 10 s eV/cm 
4.1 x 10 s eV/cm 

4.41 A 
98 meV 

48 meV 
50 meV 



TABLE I: Material parameters for single- layer M0S2. Un- 
less otherwise stated, the parameters have been calculated 
from first-principles as described in the Sees. II and III. The 
r/K-subscript on the optical deformation potentials, indicate 
couplings to the intra/intervalley phonons. 
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FIG. 6: (Color online) Scattering rates for the different 
electron-phonon couplings as a function of carrier energy at 
T — 300 K. The scattering rates have been calculated from 
the first-principles electron-phonon couplings in Figs. [3] and [S] 
The (black) dashed line shows the scattering rate obtained 
using the fitted deformation potential parameters defined in 
Eqs. ([6]) and The kinks in the curves for the optical scat- 
tering rates mark the onset of optical phonon emissions. 



for the anisotropy and the full q-dependence of the first- 
principles electron-phonon couplings. However, as mo- 
mentum and energy conservation limit phonon scattering 
to involve phonons in the vicinity of the r/i^-point, 33 
the fitting procedure yields deformation potentials close 
to the direction averaged q — > T/K limiting behavior 
of the first-principles coupling matrix elements. For the 
zero-order deformation potentials, the sampling of the 
coupling matrix elements away from the T-point in the 
k' sum in Eq. (|5P|) . leads to a deformation potentials 
which is slightly smaller than the T-point value of the 
coupling matrix element. 

The total scattering rate resulting from the fitted de- 
formation potentials is shown in Fig. [5] to give an al- 
most excellent description of the first-principles scatter- 
ing rate. With the mobility in the relaxation time ap- 
proximation given by the energy- weighted average of the 
relaxation time (see Eq. (fT5]0. the difference between the 
associated mobilities will be negligible. Hence, the defor- 
mation potentials provide well-founded electron-phonon 
coupling parameters for low-field studies of the mobility. 



C. Mobility 



general only of minor importance for the low-field mo- 
bility. This is also the case here, where it is an order 
of magnitude smaller than the other scattering rates for 
almost the entire plotted energy range. The jumps in 
the curves for the optical scattering rates at the optical 
phonon energies haj\ are associated with the threshold 
for optical phonon emission where the carriers have suf- 
ficient energy to emit optical phonons 



B. Deformation potentials 

In this section, we determine the deformation potential 
parameters in single-layer M0S2. They can be used in 
the following study of the low-field mobility within the 
Boltzmann approach outlined in Section IIVI 

The energy dependence of the first-principles based 
scattering rates in Fig. [6] to a high degree resembles that 
of the analytic expressions for the deformation poten- 
tial scattering rates in Eqs. (|2~51) . (f5T|) and (j52")l . For 
example, the acoustic and zero-order deformation po- 
tential scattering rates are almost constant in the plot- 
ted energy range. The first-principles electron-phonon 
couplings can therefore to a good approximation be de- 
scribed by the simpler isotropic deformation potentials 
in Eqs. (j6)) and (JT)). The deformation potentials are ob- 
tained by fitting the associated scattering rates for each 
of the intra and intervalley phonons separately to the 
first-principles scattering rates. The resulting deforma- 
tion potential values are summarized in Tab.|TJ In analogy 
with deformation potentials extracted from experimental 
mobilities, the theoretical deformation potentials repre- 
sent effective coupling parameters that implicitly account 



The high density of states in the K, if'-valleys of the 
conduction band in general results in non-degenerate car- 
rier distributions in single-layer M0S2. This is illustrated 
in the left plot of Fig. [7] which shows the carrier density 
versus the position of the Fermi level for different temper- 
atures. At room temperature, carrier densities in excess 
of ~ 8 x 10 12 cm -2 are needed to introduce the Fermi 
level into the conduction band and probe the Fermi-Dirac 
statistics of the carriers. Thus, only at the highest re- 
ported carrier densities of n ~ 10 13 cm -2 ^ are the carri- 
ers degenerate at room temperature. For the lowest tem- 
perature T — 100 K, the transition to degenerate carriers 
occurs at a carrier density of ~ 2 x 10 12 cm . The tran- 
sition from non-degenerate to degenerate distributions is 
also illustrated by the discrepancy between the full and 
dashed lines which shows the Fermi level obtained with 
Boltzmann statistics. 

The right plot of Fig. [7] shows the phonon-limited drift 
mobility calculated with the full collision integral as a 
function of carrier density for the same set of tempera- 
tures. The dashed lines represent the results obtained 
with Boltzmann statistics and the relaxation time ap- 
proximation using the expressions in Section IIVI for opti- 
cal phonon scattering and Matthiessen's rule for the to- 
tal relaxation time. The strong drop in the mobility from 
~ 2450 cm 2 V" 1 s" 1 at T = 100 K to - 400 cm 2 V" 1 s" 1 
at T — 300 K, is a consequence of the increased phonon 
scattering at higher temperatures due to larger phonon 
population of, in particular, optical phonons. The rela- 
tively low intrinsic room-temperature mobility of single- 
layer M0S2 can be attributed to both the significant 
phonon scattering and the large effective mass of 0.48 
m e in the conduction band. While the mobility decreases 
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FIG. 7: (Color online) Left: Carrier density as a function 
of the Fermi energy for different temperatures. The Fermi 
energy is measured relative to the conduction band edge E c 
Right: Mobility as a function of carrier density for the same 
temperatures. In both plots, the (black) dashed lines show the 
results obtained with Boltzmann statistics and the relaxation 
time approximation (only right plot). 



strongly with increasing temperature, it is relatively inde- 
pendent on the carrier density. The weak density depen- 
dence of the mobility originates from the energy- weighted 
average in Eq. (|20l) where the derivative of the Fermi- 
Dirac distribution changes slowly as a function of carrier 
density for non-degenerate carriers. As the Fermi energy 
is introduced into the band with increasing carrier den- 
sity, the derivative of the Fermi-Dirac distribution to a 
larger extent probes the scattering rate at higher energies 
leading to a decrease of the mobility. This effect is most 
prominent at T = 100 K where the level of degeneracy is 
larger compared to higher temperatures. 

Surprisingly, the relaxation time approximation is seen 
to work extremely well. The deviation from the full treat- 
ment at high carrier densities stems from the assumption 
of non-degenerate carriers and not the relaxation time 
approximation. The reason for the good performance 
of the relaxation time approximation shall be found in 
the in-scattering part of the collision integral in the full 
treatment. As the in-scattering part of the collision in- 
tegral vanishes for zero-order deformation potential cou- 
pling and is small compared to the out-scattering part 
otherwise, the full collision integral does not differ signif- 
icantly from the corresponding scattering rate as defined 
in Eq. (f30|) of Section UVl i.e. with the in-scattering part 
neglected. The good performance of the relaxation time 
approximation therefore seems to be of general validity 
for the phonon collision integral even in the presence in- 
clastic scattering on optical phonons. 

Finally, we study the temperature dependence of the 
mobility in more detail. In general, room temperature 
mobilities are to a large extent dominated by optical 
phonon scattering. This is manifested in the temperature 
dependence of the mobility which follows a \x ~ T^ 1 law 
where the exponent 7 depends on the dominating scat- 
tering mechanism. For acoustic phonon scattering above 
the Bloch-Griineisen temperature, the temperature de- 
pendence of the scattering rate in Eq. (j2l))) results in 
7 = 1. At higher temperatures where optical phonon 
scattering starts to dominate, the mobility acquires a 
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FIG. 8: (Color online) Mobility vs temperature. For compari- 
son, the mobility in the presence of only acoustic deformation 
potential scattering with the fx ~ T _1 temperature depen- 
dence is shown. The (gray) shaded area shows the variation 
in the mobility associated with a 10% uncertainty in the cal- 
culated deformation potentials. 



stronger temperature dependence with 7 > 1. In this 
regime, the exponent depends on the optical phonon fre- 
quencies and the electron-phonon coupling strength. In 
an early study of the in-plane mobility of bulk M0S2, 14 
the measured room-temperature exponent of 7 ~ 2.6 was 
found to be consistent with scattering on the homopolar 
mode via a zero-order deformation potential. 

In Fig. [51 we show the temperature dependence of the 
mobility at n = 10 11 cm~ 2 calculated with the full col- 
lision integral. For comparison, the mobility limited by 
acoustic phonon scattering with the fi ~ T^ 1 tempera- 
ture dependence is also included. To illustrate the effect 
of an uncertainty in the calculated deformation poten- 
tials, the shaded area shows the variation of the mo- 
bility with a change in the deformation potentials of 
±10%. The calculated room-temperature mobility of 
~ 410 cm 2 V -1 s _1 is in fair agreement with the recently 
reported experimental value of ~ 200 cm 2 V -1 s _1 in 
the top-gated sample of RefH where additional scattering 
mechanisms as e.g. impurity and surface-optical phonon 
scattering must be expected to exist. Over the plotted 
temperature range, the mobility undergoes a transition 
from being dominated by acoustic phonon scattering at 
T = 100 K to being dominated by optical phonon scat- 
tering at higher temperatures with a characteristic ex- 
ponent of 7 > 1. At room temperature, the mobility 
follows a temperature dependence with 7 = 1.69. The 
increase in the exponent is almost exclusively due to the 
optical zero-order deformation potential couplings and 
the Frohlich interaction while the first-order deformation 
potential couplings contribute only marginally. The ex- 
ponent found here is considerably lower than the above- 
mentioned exponent of 7 ~ 2.6 for bulk M0S2, indicat- 
ing that the electron-phonon coupling in bulk and single- 
layer M0S2 differ. Indeed, the transition from an indirect 
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band gap in bulk M0S2 to a direct band gap in single- 
layer M0S2 which shifts the bottom of the conduction 
band from the valley located along the T-K path to the 
K, if -valleys, could result in a change in the electron- 
phonon coupling. 

In top-gated samples as the one studied in Ref|H, the 
sandwiched structure with the M0S2 layer located be- 
tween substrate and gate dielectric, is likely to result in a 
quenching of the characteristic homopolar mode which is 
polarized in the direction normal to the layer. To address 
the consequence of such a quenching, the mobility in the 
absence of the zero-order deformation potential originat- 
ing from the coupling to homopolar mode is also shown 
in Fig. [SI Here, the curve for the quenched case shows a 
decrease in the characteristic exponent to 7 = 1.52 and 
an increase in the mobility of ~ 70 cm 2 V -1 s _1 at room 
temperature. Despite the significant deformation poten- 
tial of the homopolar mode, the effect of the quenching 
on the mobility is minor. 



VI. DISCUSSIONS AND CONCLUSIONS 

Based on our finding for the phonon-limited mobility in 
single-layer M0S2, it seems likely that the low experimen- 
tal mobilities of ~ 1 cm 2 V -1 s -1 reported in Refsj^HH 
are dominated by other scattering mechanisms as e.g. 
charged impurity scattering. The main reason for the in- 
crease in the mobility to 200 cm 2 V -1 s" 1 observed when 
depositing a high-K dielectric on top of the M0S2 layer 
in ReflU is therefore likely to be impurity screening. The 
quenching of the out-of-plane homopolar phonon only led 
to a minor increase in the mobility and cannot alone ex- 
plain the above-mentioned increase in the mobility. It 
may, however, also contribute. A comparison of the tem- 
perature dependence of the mobility in samples with and 
without the top-gate structure can help to clarify the ex- 
tent to which a quenching of optical phonons contributes 
to the experimentally observed mobility increase. With 
our theoretically predicted room-temperature mobility of 
410 cm 2 V" 1 s _1 , the observed enhancement in the mo- 
bility to 200 cm 2 V -1 s _1 induced by the high-ft dielec- 
tric, suggests that dielectric engineering^ is an effective 
route towards phonon-limited mobilities in 2D materials 
via efficient screening of charge impurities. 

A rough estimate of the impurity concentrations re- 
quired to dominate phonon scattering can be inferred 
from the phonon scattering rate of r _1 ~ 10 13 s _1 in 
Fig. |51 The scattering rate is related to the mean-free 
path A of the carriers via A = t(v) where (v) is their mean 
velocity. Using the velocity of the mean energy carriers 
for a non-degenerate distribution where (£k) = k^,T ', we 
find a mean-free path of A ~ 14 nm at T = 300 K. In 
order for impurity scattering to dominate, the impurity 
spacing must be on the order of the phonon mean-free 
path or smaller. This results in a minimum impurity 
concentration of ~ 5 x 10 11 cm~ 2 in order to dominate 
phonon scattering. The high value of the estimated im- 



purity concentration needed to dominate phonon scat- 
tering is in agreement with the experimental observation 
that low-mobility single-layer M0S2 samples are heavily 
doped semiconductors. 5 

As discussed in Section UH independent GW quasi- 
particle calculation a 53 ' 54 suggest that the ordering of the 
valleys at the bottom of the conduction band might not 
be as clear as predicted by DFT. If this is the case, the 
satellite valleys inside the Brillouin zone must also be 
taken into account in the solution of Boltzmann equa- 
tion. This gives rise to additional intervalley scattering 
channels that together with the larger average effective 
electron mass m* ~ 0.6 m e of the satellite valley will 
result in mobilities below the values predicted here. 

In conclusion, we have used a first-principles based 
approach to establish the strength and nature of the 
electron-phonon coupling and calculate the intrinsic 
phonon-limited mobility in single-layer M0S2. The cal- 
culated room-temperature mobility of 410 cm 2 V -1 s _1 
is to a large extent dominated by optical deformation 
potential scattering on the intravalley homopolar and 
intervalley LO phonons as well as polar optical scat- 
tering on the intravalley LO phonon via the Frohlich 
interaction. The mobility follows a fj, ~ y-i-69 t em _ 
perature dependence at room temperature characteris- 
tic of optical phonon scattering. A quenching of the 
homopolar mode likely to occur in top-gated samples, 
results in a change of the exponent in the tempera- 
ture dependence of the mobility to 1.52. With effective 
masses, phonons and measured mobilities in other semi- 
conducting metal dichalcogenides being similar to those 
of M0S2 ) 14 ' 23 room-temperature mobilities of the same 
order of magnitude can be expected in their single-layer 
forms. 
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Appendix A: First-principles calculation of the 
electron-phonon interaction 

The first-principles scheme for the calculation of the 
electron-phonon interaction used in this work is outlined 
in the following. Contrary to other first-principles ap- 
proaches for the calculation of the electron-phonon inter- 
action which are based on pseudo potentials ^ r 43 ' 55 the 
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present approach is based on the PAW methoo=, and is 
implemented in the GPAW DFT packaged— 

Within the adiabatic approximation for the electron- 
phonon interaction, the electron-phonon coupling matrix 
elements for the phonon qA and Bloch state k can be 
expressed as 



ffkq 



2MNu, 



qA 



V2V(r)|k>, (Al) 



where N is the number of unit cells, M is an appropri- 
ately defined effective mass, w qA is the phonon frequency, 
and the q-dependent derivative of the effective potential 
for a given atom a is defined by 



V^(r) = ]Te 



iq-Ri 



V Q ,V(r), 



(A2) 



where the sum is over unit cells in the system and V a iV is 
the gradient of the potential with respect to the position 
of atom a in cell R/ . The polarization vectors e qA of the 
phonons appearing in Eq. (|A1|) are normalized according 
to£a(M*/M)|eg*| a = l. 

In order to evaluate the matrix element in Eq. (|A1I) . 
the Bloch states are expanded in LCAO basis orbitals 
|iR„) where i = (a, fi) is a composite index for atomic 
site and orbital index and R„ is the lattice vector to the 
n'th unit cell. In the LCAO basis, the Bloch states are 
expanded as |k) = ^c^ik) where 



i ik ) = ^E e,;k ' R "i iR « 



(A3) 



are Bloch sums of the localized orbitals. Inserting the 
Bloch sum expansion of the Bloch states in the matrix 
element in Eq. (|Alj) . the matrix element can be written 



<k + q|e qA .V q F(r)|k) = ]T c* Cj (tk + q|e qA • V q V(r)|jk> 

i(k+G)-(R„ R m )-»q- (R m — R; ) 

(A4) 



=^E c ^E e " 

ij hnn 

x («R m |e qA • V;V(r)|jR„) 



where the k-labels on the expansion coefficients have 
been discarded for brevity. The phase factor from the re- 
ciprocal lattice vector G can be neglected since R n ■ G = 
2itN. By exploiting the translational invariance of the 
crystal the matrix elements can be obtained from the 
gradient in the primitive cell as 



(iR m |e qA - V;F(r)|jR n ) 

= (iRm - R,|e qA ■ V Q V{r)\jR n 



R) 



(A5) 



by performing the change of variables r' = r — R; in 
the integration. Inserting in Eq. (|A4[) and changing the 
summing variables to R m — R; and R„ — R; we find for 
the matrix element 



k + q|e qA • VqF(r)|k) = ^_ £ c * Cj £ elk .(R„-R m )- i : q .R m (lRm |e qA . V V(v)\jR n ) 



jf E c *^ E e lk ' (R "- R - ) - tq - R -( J R m |e qA • V V(r)|jR„) 



(A6) 



Here, the sum over k in the first equality produces a 
factor of N. The result for the matrix element in Eq. (IA6|) 
is similar to the matrix element reported in Eq. (22) of 
ReflU where a Wannier basis was used instead of the 
LCAO basis. 



are obtained as 



From the last equality in Eq. (|A6I) . the procedure for 
a supercell-based evaluation of the matrix element has 
emerged. The matrix element in the last equality in- 
volves the gradients of the effective potential VoV^ with 
respect to atomic displacement in the reference cell at Ro . 
These can be obtained using a finite-difference approxi- 
mation for the gradient where the individual components 



dV 

dR a i 



V + (r) - V~(r) 
28 



(A7) 



Here, denotes the effective potential with atom a 
located at R Q in the reference cell displaced by ±5 in 
direction i = x,y, z. The calculation of the gradient thus 
amounts to carrying out self-consistent calculations for 
six displacements of each atom in the primitive unit cell. 
Having obtained the gradients of the effective potential, 
the matrix elements in the LCAO basis of the super- 
cell must be calculated and the sums over unit cells and 
atomic orbitals in Eq. (IA6[) can be evaluated. 

In general, the matrix elements of the electron-phonon 
interaction must be converged with respect to the super- 
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cell size and the LCAO basis. In particular, since the 
supercell approach relies on the gradient of the effective 
potential going to zero at the supercell boundaries, the 
supercell must be chosen large enough that the potential 
at the boundaries is negligible. For polar materials where 
the coupling to the polar LO phonon is long-ranged in 
nature, a correct description of the coupling can only be 
obtained for phonon wave vectors corresponding to wave- 
lengths smaller than the supercell size. Large supercells 
are therefore required to capture the long-wavelength 
limit of the coupling constant (see e.g. main text). 



1. PAW details 

In the PAW method, the effective single-particle DFT 
Hamiltonian is given by 

H=-lw 2 +v cS (r)+Y / T, \PV A KiM 2 l (A8) 



where the first term is the kinetic energy, v e e is the effec- 
tive potential containing contributions from the atomic 
potentials and the Hartree and exchange-correlation po- 
tentials and the last term is the non-local part includ- 
ing the atom-centered projector functions pf(r) and 
the atomic coefficients AH"^. In contrast to pseudo- 
potential methods where the atomic coefficients are con- 
stants, they depend on the density and thereby also on 
the atomic positions in the PAW method. 

The diagonal matrix elements (i.e. q = 0) in Eq. (|A4|) 
correspond to the first-order variation in the band en- 
ergy £k upon an atomic displacement in the normal mode 
direction e q >. Together with the matrix elements for 
q ^ 0, these can be obtained from the gradient of the 
PAW Hamiltonian with respect to atomic displacements. 
The derivative of (IA8|) with respect to a displacement of 
atom a in direction i results in the following four terms 



dH 

8R a i 



dr 



dR 



dRai 



dR a i 
r)/\ TT a 

\p a ) *lfe/p? 

lP%1 ' dR m ^* 



(A9) 



While the gradient of the projector functions in the 
first and second terms inside the square brackets can 
be evaluated analytically, the gradients of the effective 
potential and the projector coefficients in the last term 
are obtained using the finite-difference approximation in 
Eq. (1A71) . 

Once the gradient of the PAW Hamiltonian has been 
obtained, the matrix elements in Eq. (|A4[) can be ob- 
tained. Under the assumption that the smooth pseudo 
Bloch wave functions "0k from the PAW formalism is a 



good approximation to the true quasi particle wave func- 
tion, the matrix element follows as 



(k + q|e qA 



VqV(r)|k) 



(V'k+qleqA ' V q i7|^k 



(A10) 



where V q i? is given in Eq. (IA2|) . With the pseudo wave 
function expanded in an LCAO basis, the matrix element 
is evaluated following Eq. (|A6|) . 

In order to verify the calculated matrix elements, we 
have carried out self-consistent calculations of the vari- 
ation in the band energies with respect to atomic dis- 
placements along the phonon normal mode. For the 
coupling to the T-point homopolar phonon in single- 
layer M0S2, we find that the calculated matrix element 
of 4.8 eV/A agrees with the self-consistently calculated 
value to within 0.1 eV/A. 

As the matrix elements of the electron-phonon inter- 
action are sensitive to the behavior of the potential and 
wave functions in the vicinity of the atomic cores, vari- 
ations in electron-phonon couplings obtained with dif- 
ferent pseudo-potential approximations can be expected. 
We have confirmed this via self-consistent calculations of 
the band energy variations, which showed that the cou- 
pling to the homopolar L-point phonon increases with 
0.3 eV/A compared to the PAW value when using norm- 
conserving HGH pseudo potentials. 



Appendix B: Frohlich interaction in 2D materials 

In the long-wavelength limit, the lattice vibration of 
the polar LO mode gives rise to a macroscopic polar- 
ization that couples to the charge carriers. In three- 
dimensional bulk systems, the coupling strength is given 
by Eq. [5] which diverges as l/q. Using dielectric contin- 
uum models, the coupling has been studied for confined 
carriers and LO phonons in semiconductor heterostruc- 
ture o 44 ' 45 where a non-diverging interaction is found. For 
atomically thin materials, however, macroscopic dielec- 
tric models are inappropriate. Using an approach based 
on the atomic Born effective charges Z* , a functional 
form that fits the calculated values of the Frohlich inter- 
action in Fig. [3] is here derived. 



1. Polarization field and potential 

In two-dimensional materials, the polarization from 
the lattice vibration of the polar LO phonon is oriented 
along the plane of the layer. It can be expressed in terms 
of the relative displacement u q of the unit cell atoms as 



Pq(*) = 



"T u qJq 



/q(*)> 



(Bl) 



where q is the two-dimensional phonon wave vector, 
is the optical dielectric constant, Z* is the Born effective 
charge of the atoms (here assumed to be the same for all 
atoms) and / q describes the profile of the polarization in 
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the direction perpendicular (here the z-direction) to the 
layer. The associated polarization charge p = — V • P 
is given by p q = — iq ■ P q . The resulting scalar poten- 
tial which couples to the carriers follows from Poisson's 
equation. Fourier transforming in all three directions, 
Poisson's equation takes the form 

Z* 



(q 2 + k 2 )^(k) 



£00 A 



q • Vq( fc ) 



(B2) 



where k is the Fourier variable in the direction perpen- 
dicular to the plane of the layer and q || u q for the LO 
phonon. 

In 2D materials, the z-profile of the polarization field 
can to a good approximation be described by a 5- 
function, i.e. 



(B3) 



Inserting the Fourier transform / q (fc) = 1 in Eq. (|B2[) . 
we find for the potential 



*q(*) 



,Z*U q 

.Z*u. 



dk e* fcz / q (*0 — 



k 2 



£<->o A 



R e -8l*l 



(B4) 



which in agreement with the findings of Refs l44ll45l does 
not diverge in the long- wavelength limit. 



Appendix C: Evaluation of the inelastic collision 
integral 

Following Eq. (|2"o]) in the main text, the inelastic colli- 
sion integral for scattering on optical phonons is split up 
in separate out- and in-scattering contributions, 



l'm e l[<Pk] = I'inel{<f>k] + I'^eAfa] 



(CI) 



which are given by Eq. (|27|) and (|29[) , respectively. With 
the assumption of dispersionless optical phonons, the 
Fermi-Dirac and Bose-Einstein distribution functions do 
not depend on the phonon wave vector and can thus be 
taken outside the q-sum. 

The out-scattering part of the collision integral then 
takes the form 



I'° 



2tt 



■>k\ = 



^a°(1 - /k +q ) E l^l 2 *( £ k+q - e k - 



+ (l + 7V A )(l-/ k °_ q )]>> q | 2 % k 
q 



£k 



(C2) 



2. Electron-phonon interaction 



For the in-scattering part, the additional cos$k,k±q fac- 
tor is rewritten as 



In three-dimensional bulk systems, the g-dependence 
of the Frohlich interaction is given entirely by the 1 jq- 
divergence of the potential associated with the lattice 
vibration of the polar LO phonon. However, in two di- 
mensions, the interaction follows by integrating the po- 
tential with the square of the envelope function x(z) of 
the electronic Bloch state. Hence, the Frohlich interac- 
tion is given by the matrix element 



Slo(<?) = dz x*(z)(f>q(z)x(z). 



(B5) 



For simplicity, we here assume a Gaussian profile for the 
envelope function 

1 



xO) 



rl/4 



-*72ff a 



(B6) 



where a denotes the effective width of the electronic 
Bloch function. With this approximation for the enve- 
lope function, the Frohlich interaction becomes 



5Lo(<z)=.9FrXe« 2ff2 / 4 erfc( g( 7/2) 
« 3Fr x erfc(gcr/2) 



(B7) 



where $Fr is the Frohlich coupling constant, erfc is the 
complementary error function and the last equality holds 
in the long- wavelength limit where q^ 1 a. As shown 
in Fig. [5] in the main text, this functional form for the 
Frohlich interaction gives a perfect fit to the calculated 
electron-phonon coupling for the polar LO mode. 



cos fl k ,k±q 



k ± gcos #k,< 



kjl± 



(C3) 



leading to the following general form of the in-scattering 
part of the collision integral 



incl L 



_ 2tt 1 

^(Wk+q)- 



<t>k+q 



kjl 



x 2J l5q| 2 (k + gcos6»k, q )<5(£k+q - £k - f^x) 
q 



+ (l + A A )(l-/ k °_ q )- 



kjl 



X V* |5q| 2 ( fc - <7COs6> k ,q)<5(£k-q ~ £k + fuV\) 



q 



(C4) 



This depends on the deviation function in <pk±q = </>(£k± 
and thus couples the deviation function at the initial 
energy £k of the carrier to that at £k ± Swa- The two 
terms inside the square brackets account for emission and 
absorption out of the state k ± q and into the state k, 
respectively. 
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1. Integration over q 

The in- and out-scattering contributions to the inelas- 
tic part of the collision integral in Eqs. (|C2I) and (IC4[) 
can be written on the general form 

^ OUt (k) = J2 /(?)*(ek±q - ^ T (C5) 



where the function / accounts for the q-dependent func- 
tions inside the q-sum. The conservation of energy and 
momentum in a scattering event is secured by the 5- 
function entering the collision integral. 

Following the procedure outlined in Refill, the integral 
over q resulting from the q-sum in the collision integral is 
evaluated using the following property of the (^-function, 



[ dq F{ q )5[g{ q )]=Y.Tl 



F(q n ) 



(C6) 



Here, F(q) = qf{q) and q n are the roots of g. Rewriting 
the argument of the (^-function in Eq. (|C5I) as 



h 2 (k±q) 2 h 2 k 2 _ 

£k± q ~£kT huq X = — — - T ^qA 



2m* 



h 2 

2m* 

h 2 

2m* 



q 2 ± 2fcgcos#k,q T 



9(q) 



we get 



J'(k) = 



.4 



{2-Kf 



dO 



k.q 



2TO* X - gn/(gn) 

V ^ W(qn)\ 



(C7) 



(C8) 



Depending on the q-dependence of the phonon frequency, 
different solutions for the roots q n result. 



In the case of phonon emission, both possible roots in 
Ea. (IU9|) 



q±=k COS 6*k,q ± k^/cOS 2 0k, q ~ — — 



(C12) 



can take on positive values. However, the minus 
sign inside the square root restricts the allowed val- 
ues of the integration angle to the range 9^ q < 9 = 
arccos (HuJx/sk) 1 / 2 . Furthermore, in order to secure a 
positive argument to the square root, the electron en- 
ergy must be larger than the phonon energy, £k > huj\. 
Hence, the emission terms can be obtained as 



.4 



Z' e (k) =7^%e( ek - hw x ) 
(2tt) 2 



^^ifm±£m. ( ci3) 



cos 2 9k,, 



2. Collision integral for optical phonon scattering 

In the following two sections, analytic expressions for 
collision integral in the case of zero- and first-order op- 
tical deformation potential coupling are given. For the 
Frohlich interaction the integration over 9 can not be re- 
duced to a simple analytic form and is therefore evaluated 
numerically. 



a. Dispersionless optical phonons 



a. Zero-order deformation potential 



With the assumption of dispersionless optical phonons, 
i.e. LUq\ — lu\, the roots of g in Eq. (|C7[) become 



q/k = =FCOS#k,q± J cos 2 #k, q ± — ^ (C9) 
V £k 

for absorption (upper sign in the first and last term) and 
emission (lower sign in the first and last term), respec- 
tively. 

For absorption there is only one positive root q°L given 
by the plus sign in front of the square root, 



q a + = -/scos(9 k q + £;Wcos 2 6>kq + — (CIO) 

V £k 

where all values of #k,q are allowed. The absorption terms 
in the collision integral then becomes 



m 

q kh 2 



COS 2 # k ,q + ^ 



(Gil) 



For the zero-order deformation potential coupling in 
Eq. (J7J, the ^-integration in the out-scattering part of 
the collision integral yields a factor of 2tt resulting in the 
following form 



wOllt , 

1 inellWJ 



m*Dj x <j) k 
2h 2 puxl-f° 



(1 - /k+q) + (1 - /k -q)e^ /fcBT e( £k - HU X ) 



~(C14) 



For the in-scattering part, the (/-independence of the de- 
formation potential results in a vanishing q-sum and the 
in-scattering contribution is zero, i.e. 



(C15) 
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b. First-order deformation potential 



and 



For first-order deformation potentials, we find for the 
out-scattering part 



wOllt r , 



m* 2 ^ 2 x (j) k 



(1-/° )(2e k + ^ A ) 



+ (1 - /k- q )(2£k - fiw A )e ;iaJA/fcBT e(e k - fi^ A ) 



(C16) 

For the in-scattering part, the result for the first and 
second term inside the square brackets in Eq. (|C4[) is 



*2~2 1 _ fO 



k+q 



*2^2 1 _ fO 
H\ 1 Jk- 



cl[<Pk\ 



k-q 



n 4 pWA l-/k A _ ^ 



x (i + ivj) ( £k - fiw A ) e( Ek - /iw A ) 



1 / 7T / ^CJa / ^ A \ 

x — — + arccos \/ harcsinW , (C18) 

7T \ 2 V £k V £k I 



(C17) which accounts for absorption and emission, respectively. 
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